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Abstract 

In modeling strong gravitational lens systems, one often adopts simple models, such as singular isother- 
mal elliptical density plus lowest-order external perturbation. However, such simple models may mislead us 
if the real mass distribution is more complicated than that in the assumed models. In particular, assump- 
tions on mass models are crucial in studying flux ratio anomalies that have been suggested as evidence for 
a cold dark-matter substructure. We reinvestigated four quadruple lens systems using power-law Fourier 
models, which have advantages of clear physical meanings and the applicability of a linear method, as well 
as a simple singular isothermal elliptical density model. We also investigated the effect of external per- 
turbations, including a singular isothermal sphere, lowest-order expansion, and next-order expansion. We 
have found that the cos 36 terms of the primary galaxy and/or of external perturbation significantly reduce 
the x 2 in PG 1115+080 and B1422+231. In particular, we could reproduce the flux ratios of B1422+231 
with next-order external perturbation assuming 5% flux uncertainties, suggesting that external perturba- 
tion cannot be described by a simple singular isothermal sphere approximation. On the other hand, we 
could not fit B0712+472 and B2045+265 very well even with our models, although the \ 2 were reduced 
, compared with the case of using simple models. Our results clearly demonstrate that both the primary 

lens galaxy and the external perturbation are often more complicated than we usually assume. 
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Q ■ 1. Introduction 

Gravitational lenses of distant quasars are quite useful in cosmological and astrophysical studies, such as measuring 
the Hubble constant from time delays (Refsdal 1964), studies of galaxy formation and evolution, and direct probes of 

. galaxy structure, including cold dark-matter substructures. The method used to measure the Hubble constant does 
k> ' n °t re ly 011 the distance ladder, which inevitably suffers various systematic errors. After the discovery of the first 
lens system, Q 0957+561 (Walsh et al. 1979), many lens systems have been studied in detail, adopting various lens 

■ models. However, the model uncertainties, including several intrinsic degeneracies (Falco et al. 1985; Gorenstein et 
al. 1988; Saha 2000), prevent us from accurately deriving the Hubble constant (e.g., Wambsganss, Paczyhski 1994; 
Keeton et al. 1997; Kccton, Kochanck 1997; Courbin et al. 1997; Wucknitz, Refsdal 1999; Witt et al. 2000; Chiba, 
Takahashi 2002; Wucknitz 2002; Oguri et al. 2002; Oguri, Kawano 2003). 

Among several applications of gravitational lenses, flux-ratio anomalies have attracted much attention. There are 
mainly three explanations for them: Effects of the interstellar medium, problems in smooth models, and lensing by the 
substructure (e.g., Kochanek, Dalai 2003). While substructure lensing is often regarded as being the best explanation 
for the flux anomalies (e.g., Mao, Schneider 1998; Chiba 2002), the former two explanations have not yet been 
sufficiently investigated to reach the conclusion that the flux-ratio anomalies are connected with the substructure. 
The limitation comes from the fact that it is difficult to apply all of the models, which are physical and complicated, 
to the observed lens systems. In addition, since galaxies arc believed to consist of not only luminous matter, but 
also dark matter, the shapes of lens galaxies do not have to be elliptical, nor other simple shapes, as the isophotcs 
indicate. Therefore, it is of great importance to investigate the contribution of many forms of the angular structure 
of the lens potential to find the true origin of the flux anomalies. Indeed, Evans and Witt (2003) showed that some of 
lens systems that have been claimed to be lenses with flux anomalies can be well fitted by a smooth lens model using 
a power-law Fourier potential (see also Oguri, Kawano 2003). Kochanek and Dalai (2003) studied the contribution of 
only the isothermal m = 3 and m = 4 multipolcs to the flux ratio anomalies reported so far. We then reinvestigate 
some of the systems with not only multipolcs, but also several types of external perturbations. 

Our aim in this paper is to discuss the contribution of both multipoles of the primary lens galaxy and external 
perturbation of the group, in particular the cos 3$ terms that are not usually included in lens modeling. Wc consider 
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whether some multipoles can produce an acceptable fit without substructure or not. We also discuss the quadruple 
lens system B0712+472 in detail because of the interesting observations that the ellipticity of the galaxy and the shear 
of the group seem to be significantly misaligned (~ 90°) to each other, and that previous authors estimated the shear 
to be small compared to observations using the one-shear model (Jackson et al. 1998; Fassnacht, Lubin 2002). Since 
more than three images are necessary to constrain the mass models strongly, we restrict our considerations to four 
quadruple lens systems (PG 1115+080, B1422+231, B0712+472, and B2045+265). There are three reasons for this 
restriction. First, the existence of a substructure of the lens systems has been claimed to account for the observed flux 
ratios (e.g., Kochanek, Dalai 2003). Second, there are no satellite galaxies that complicate the lens model. Third, there 
are sufficient astrometric, photometric, and spectroscopic (redshift) data. In this paper, we consider only isothermal 
models; however our analysis can be easily extended to non-isothermal radial mass profiles. 

This paper is organized as follows. In section 2, we describe the models that we used to investigate four systems. 
In section 3, we show how we applied the models to four systems and give the results. We finally discuss the results 
and give conclusions in section 4. The linear method that we used is presented in Appendix 1. 

2. Lens Model and Method 

2.1. Fundamental Equations of Gravitational Lensing 

A source at x s and an image i at Xi are related through the lens equation, 

x t - x s = Vil>(xi), (1) 
where tp(x) is the projected lens potential. Thus, the magnification factor is given by 

dx, 
dx 

We can also calculate the equation of the relative time delay between image i and j using 
1 + z d D d D s 



li{xi) = 



dct 



2) 



Atij = 



c D ds 

\Ax 2 --Ax 2 -^{x. i ) + ^{x j ) 



(3) 



where is the redshift of the lens object and c is the speed of light; D^, D s , D^s are the angular diameter distances 
to the lens, to the source, and from the lens to the source, respectively; Ax, = x- t — x s (Schneider et al. 1992). The 
two-dimensional Poisson equation relates the dimensionlcss density distribution k(x) to the lens potential as 

VV(a;) = 2k(x). (4) 
2.2. Lens Model 

For the potential of the primary lens galaxy, we apply the singular isothermal elliptical density (SIED) model and/or 
the power-law Fourier potential (PLF; Oguri, Kawano 2003; Evans, Witt 2003) model. The former model may be 
called the "standard" model in the sense that it is widely used in the mass modeling of a lens galaxy. The two- 
dimensional density profile of SIED is given, in the Cartesian coordinates (x, y) centered on the primary lens galaxy, 
as 

K^y{x)=b{x 2 +y 2 lq 2 )- 1 l\ (5) 

where b is the mass parameter that characterizes the Einstein radius and q is the axis ratio. The ellipticity, e, is 
given by 1 — q. We used the latter model, PLF, to investigate the contribution of multipoles, particularly an m = 3 
multipole, which is neglected in usual lens modelings of lens systems. In addition, this model has an advantage that 
we can apply a linear method (Appendix 1). The lens potential of the power-law Fourier model is given by 

i> sa i a xy = r f3 F(6), (6) 

where a Fourier expansion of the angular function F{9) is 

oo 

F{9) = a 4- ^ (a m cosm6 + & m sinm0) (7) 

rn — l 

oo 

= a + A m coa{m(0-8 m )}, (8) 

rn — l 
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in polar coordinates (r,6). The variables A m = (a 2 n + b 2 n ) 1 ' 2 and 9 m are the amplitude and the position angle of the 
multipole, respectively. In the present work, we considered only the case of j3 = 1, which has a similar radial profile as 
in SIED, while the parameter (3 can take values between (point mass) and 2 (mass sheet). 

While the influence of the galaxy group which embeds the primary lens galaxy should be taken into account (e.g., 
Keeton et al. 1997; Witt et al. 2000), detailed information on the density distribution of the group is uncertain in most 
cases. The potential near the lens galaxies is often described by a singular isothermal sphere (SIS) or the potential 
expanded around the center of the lens galaxy (Kochanek 1991; Bernstein, Fischer 1999; Keeton et al. 2000) as 

0ex = 00 + 01 x 
i J* 



1 1 

2«o + ^7cos2( 



'1j 



-r 3 



— ersin(# — a ) — — (5sin3( 
4 6 



+0(r A ). (9) 

The first three terms cannot be fully constrained by observables, because of the lens degeneracies (Falco et al. 1985; 
Gorenstein et al. 1988; Saha 2000). Thus, one usually uses an external shear (XS), 



0xs = |r 2 cos2(6l-(9 7 ) = y r 2 cos26> + y r 2 sin 20, (10) 



where 71 = 7cos20 7 and 72 = 7sin20 7 . In the present work, we also added the terms of third order (r 3 ) to the lens 
potential in a similar way of e.g., Bernstein and Fischer (1999) enforcing S = — (3/2)er and 8s = 9 a (as for a singular 
isothermal group). In this case, the explicit expression of the additional potential (X3) is 
3 

0x3 = y<r [sfn(0 - a ) + sin3(0 - a )] . (11) 

In summary, we used SIS, SIED, and PLF models for the primary lens galaxy, and XS, SIS, or (XS+X3) models for 
the group. 

2.3. Method 

We used image positions and flux ratios as constraints on the lens models. The relative time delays between images 
were used not only for determining the Hubble constant, but also for additional constraints on the lens models. The 
number of constraints was more than ten for a quadruple lens system, while the typical number of model parameters 
is equal to or more than seven. The goodness of fit was determined by x 2 statistics based on the image positions, the 
flux ratios, and the relative time delays: 

L ' ( 12 ) 



x 2 



pos 



(model) \ 2 



xl^E^M^ 



fij 



Xtime 



^(At„--A^ lod0l) ) 2 

e - y j — -■ ( u ) 



In the above equations, Xi are the observed image positions, are the observed flux ratios, and Atij are the observed 
relative time delays. The suffix "(model)" indicates the model predictions for each observable. The uncertainties in 
the observables are represented by a X i, a fij, and a t ij- We minimized the total chi-square, Xtotai = Xp OS + x| U x + Xtime 
(Press et al. 1992), to obtain the best-fit parameters, itbost- In this fitting, (i) we solved non-linear lens equations 
to derive the model predictions for the image positions, (ii) calculated the relative time delays and the flux ratios 
from the predicted image positions and the model parameters, and (iii) evaluated and minimized Xtotai- Finally, we 
evaluated the reduced chi-square, Xtotai /-^dof, where A^jof is the number of degrees of freedom. In the case of the PLF 
model, the \ 2 minimization was significantly simplified because of the linearity of the equations (see Appendix 1). 

3. Results 



About 80 gravitational-lens systems have been reported so far. One third of them arc quadruple-lens systems, 
and the rest are double-lens systems. However, we considered only quadruple lens systems, because more than six 
observables are typically needed to strongly constrain lens models. Specifically, we analyzed the systems PG 1115+080 
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Table 1. Results of PG 1115+080. 



Model 



Best-fit parameters 



Note 



SIED 

SIS+XS 

SIED+XS 

SIS+(XS+X3) 

SIS+SIS 

(m=2)+XS 



537.9/8 
257.0/8 

(27.4+1. 5+1.8)/6 

(0.6+2.9+0.3)/6 

(2.6+2. 7+0. 2)/7 
(27.8+1. 5+1. 7)/6 



(m=3)+XS (1.0+3.3+0. 7)/6 



Notes. The x 2 is the total chi-square, x 2 a 



e = 0.32, 6 C = 66° 

7 = 0.11, 7 = -115° 

e = 0.18, 6> c = 89° 

7 = 0.09,6> 7 = -132° 

7 = 0.12,6> 7 = -115° 

a = 0.008, 6 a = -160° 

r — i o"4 — — 1 1 fi c 

'group ±J ^- u group iiu 

7 = 0.11,6> 7 = -129° 
A 2 /a = O.O3,0 2 = 6° 
7 = 0.12,6> 7 = -115° 
A 3 /a = 0.002, 6» 3 = -14° 



Acceptable fit 
Acceptable fit 

Acceptable fit 



-v? . Model 



denotes the used potential of galaxy+group. The position of the 
group center is (r grc ,up , groI1 p ) . 



(Weymann et al. 1980), B1422+231 (Patnaik et al. 1992), B0712+472 (Jackson et al. 1998), and B2045+265 (Fassnacht 
ct al. 1999). We adopted a flat lambda-dominated universe with (Qo,^a) = (0.3,0.7), where f2o is the density parameter 
of matter and £l\ is the dimensionless cosmological constant. The Hubble constant in units of 100 km s _1 Mpc -1 is 
denoted by h, as usual. 

3.1. Data 

The positions and fluxes of most of the gravitationally lensed images are relatively well constrained within 10 mas 
and about five percent uncertainties, respectively, while measuring time delays from light curves are very difficult. For 
most of the systems, we used the data of the image positions and the fluxes by the CfA/Arizona Space Telescope Lens 
Survey 1 (CASTLES; Falco ct al. 1999), and for some of the systems the data was taken from the papers cited below. 
The uncertainties of the flux ratios were slightly broadened (see below) to account for possible systematics, such as 
possible microlensing (e.g., Mao, Schneider 1998), extinction, and quasar variability during a time delay. 

The system PG 1115+080 is a QSO (z s = 1.72) gravitationally lensed by an early-type galaxy of Zd = 0.31, which 
looks very spherical (Keeton ct al. 1998; Iwamuro et al. 2000). The lens galaxy belongs to a group of galaxies. The 
QSO images consist of four components (fold-type) with the maximum angular separation of 2.32 arcsec. In order to 
constrain models, we used the CASTLES data of the image positions and the fluxes, and Barkana (1997)'s analysis of 
the time delays. The flux uncertainties were broadened to 20%. 

The second system, B1422+231, that we analyzed is known as a system with possible substructure microlensing 
(Mao, Schneider 1998; Keeton 2001; Chiba 2002; Bradac et al. 2002; Keeton et al. 2003; Kochanek, Dalai 2003). The 
elongated galaxy (zd = 0.34) associated with a poor group gravitationally lenses a QSO (z s = 3.62), and leads to a 
cusp-type four-imaged system. The ellipticity, e, of the galaxy is 0.27±0.13 (Impey et al. 1996). There are radio data 
of this system (Patnaik et al. 1992, 1999), having uncertainties of 50 //arcsec for QSO images and a few percents for 
their fluxes. Some other observational properties, such as the effects of the Lyman a (Bechtold, Yee 1995; Petry et 
al. 1998) and the variability (Yee, Bechtold 1996; Koopmans et al. 2003), have been observed. We used the data of 
Patnaik et al. (1999), assuming flux uncertainties of 5%. 

We also investigated an interesting system B0712+472 (Jackson ct al. 1998). The rcdshifts of the elongated primary 
lens galaxy and the source QSO are 0.41 and 1.34, respectively. Interestingly, the primary lens galaxy is not associated 
with a group of galaxies, but a foreground group with a mean redshifts of 0.29 are overlapped (Fassnacht, Lubin 2002). 
The light distribution of the galaxy is very elongated with an ellipticity of between 0.4 and 0.5. The positions of the 
QSO and the galaxy images that we used were from Jackson et al. (1998) and CASTLES, respectively. For the fluxes, 
radio data of Jackson et al. (1998) was used, assuming 20% uncertainties. 

The system B2045+265 was discovered by a radio survey of the Cosmic Lens All-Sky Survey (CLASS). Five radio 
components were clearly observed (Fassnacht ct al. 1999). The four components arc the lensed QSO (z s = 1.28) 
and another component is the radio core of the primary lens galaxy (zd = 0.87), which is associated with a compact 
group of galaxies. The optical spectrum of the galaxy is similar to that of nearby Sa galaxies. We used the data of 
CASTLES for the galaxy and that of Fassnacht et al. (1999) for the QSO. The flux uncertainties were assumed to be 
20% uncertainties. 
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Table 2. Results of B1422+231. 



Model 








Best-fit parameters 




Note 


SIED-r 


XS 




(2.5+12.2)/4 


e = 0.28, e = -56° 














7 = 0.16, 7 = 127° 






SIED-r 


SIS 




(1.5+9.0)/3 


e = 0.22, C = -58° 














1 group — tJ. group 


= 128° 




SIED-r 


■(XS- 


HX3) 


(0.1+3. 8)/2 


e = 0.25, c = -59° 




Acceptable fit 










7 = 0.20 , 7 = 134° 














a = 0.02, CT = 63° 






[SIED- 


h(m= 


=3)]+XS 


(0.1+6.1)/2 


e = 0.29, c = -55° 














7 = 0.17, 7 = 129° 














A 3 /a = 0.007, 3 = 


-16° 




(m=3) 


+XS 




781.8/4 


7 = 0.21, 7 = 124° 














A 3 /a = 0.01, 3 = - 


-18° 


Too large A 3 



Notes. The \ z is the total chi-square, Xpos + Xfl ux - Model denotes 
the used potential of galaxy+group. 



3.2. Implications from Models 

Table 1 gives the results of PG 1115+080 with typical models that we described in previous section. Some models 
[(m = 3)+XS, SIS+(XS+X3), and SIS+SIS] well reproduce the image positions, the flux ratios, and the relative time 
delays. While most previous authors applied a SIED+XS or a SIED+SIS, our results clearly imply that the lens 
requires the term cos 30, which might be a component of the primary lens galaxy (to = 3 multipolc), or be that of the 
external perturbation of the group (X3). A very small amplitude (A^/ao « 0.002) of the to = 3 multipole of the galaxy 
is consistent with the galaxy light distribution, which is very spherical. While some of the models including only XS as 
the external perturbation fail to reproduce the obscrvablcs, inclusion of the contribution from the group represented 
by SIS and XS+X3 always gives a better fit, i.e., an X3 perturbation is important to improve the fit. The best fit 
of the angles in the external perturbation are 7 = —115° and CT = —160°. These angles are consistent with optical 
(Impey et al. 1998) and X-ray (Grant et al. 2003) observations of the group. Therefore, we conclude that there is a 
degeneracy between an m = 3 multipole of the lens galaxy and an X3 perturbation of the group in this lens system, 
possibly because the galaxy is very spherical and the external shear is small. The mass distribution of the galaxy 
might be spherical and the ellipticity of the galaxy might be unnecessary. Schechter et al. (1997) also used the model 
of a SIS+SIS, but they did not allow the galaxy position to vary. Allowing a variable galaxy position is essential. 

The results of B1422+231 arc given in table 2. While the SIED+XS and even the SIED+SIS cannot reproduce 
the observablcs of the system, particularly the flux ratios, including the perturbation of X3 significantly reduced the 
X 2 /iVdof to less than 2 [see the (SIED+m = 3)+XS model]. The fact that the angles of the external perturbations are 
somewhat strange, |0 7 — CT | = 1 134° — 63° | = 71° , implies that the structure of the group is not so simple as SIS. Recent 
X-ray observations (Grant et al. 2003) support this interpretation. Moreover, the direction of the a perturbation in 
our modeling [equation (9)] coincides with the direction of the gradient of the group density distribution. Note that a 
part of the large misalignment of the external perturbations (7 and a) may arise from limiting the amplitudes and the 
directions of the r 3 perturbation [equation (11)]. Since the system was well-fitted with 5% flux uncertainties, it may 
be unnecessary to consider microlcnsing by substructure, or other compact structures, for this system, which supports 
the analysis of Keeton, Gaudi, and Petters (2003) using the cusp relation. 

We also show results of B0712+472 in table 3. In this system, it does not seem that fitting by the SIED+XS model 
was significantly improved by including other additional terms, such as an to = 3 multipole of the galaxy and X3 terms 
of the group. Therefore, it may be better to consider the existence of substructure(s) in this galaxy. However, we note 
that the curves of the flux ratios are variable, and perhaps further observational studies are needed (see Koopmans et 
al. 2003). Another interesting result is that both the ellipticity of the primary lens galaxy and the amplitude of the 
external shear of the foreground group are large in a iwo-shear model, such as SIED+XS. This is in marked contrast 
with the one-shear model studied so far (Keeton et al. 1998; Jackson et al. 1998; Fassnacht, Lubin 2002) in which the 
ellipticity of the light is larger than that of the mass and the external shear is very small (7 < 0.1). We believe that 
the results with ellipticity plus external shear are more realistic, since the two shears are in good agreement with that 
of observed light. The conflict between the two can be easily understood by the fact |0 C — 7 | rs 90°, where e is the 
angle of the mass ellipticity. Here, we roughly represent the potential of the SIED as the PLF to to = 2, 

V'sied «r[a + Ae cos 2(0- e )], (15) 
where A c is the amplitude of the mass ellipticity. Their contribution to the lens potential is approximately expressed 
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Table 3. Results of B0712+472. 



Model X/Ndof Best-fit parameters Note 

SIED (33.5+54.6)/6 e = 0.33, 9 C = 50° 

SIS+XS (53.2+86.8)/6 7 = 0.12, 9 1 = 49° 

SIED+XS (10.1+3.0)/4 e = 0.73, 6> c = 56° 

7 = 0.33, 9 1 = -34° 
SIED+SIS (9.3+4.9)/3 e = 0.61, 6 e = 55° 

^"group — 6. 0, $group — 152 

SIED+(XS+X3) (2.5+4.8)/2 e = 0.69, O = 56° 

7 = 0.32, 6» 7 = -43° 
£7 = 0.19, ff = 123° 

[SIED+(m=3)]+XS (3.9+0. 7)/2 e = 0.63, 6> c = 54° Too large A 3 

7 = 0.22, 6» 7 = -39° 

A 3 /a = 0.01, 9 3 = 37° 

Notes. The x l l& the total chi-squarc, Xpos+Xfl ux - Model denotes 
the used potential of galaxy+group. 



Table 4. Results of B2045+265. 



Model 






Best-fit parameters Note 


SIED+XS 




481.1/4 


e = 0.19, 9 C = -70° 


SIED+SIS 




(16.0+48.0)/3 


7 = 0.10, 9 1 = -69° 

r sis = l."5, 9 sis = -66° Not smooth 


SIED+(XS- 


I-X3) 


(2.3+15. 2)/2 


e = 0.81, 9 C = -68° 
7 = 0.66 , 7 = 21° 
cr = 0.21, 9 1 = 160° 


[SIED+(m= 


=3)]+XS 


311.4/2 


e = 0.37, 9 C = 24° 
7 = 0.22, 7 = -67° 
A 3 /a =0.03, 9 3 = -70° 



Notes. The x 1 ls the total chi-square, Xpos + Xfl ux - Model denotes 
the used potential of galaxy+group. 



(without r dependences) as A c cos2(0 - 9 C ) + (7/2) cos2(0 - 7 ) = A c cos2(0 - 6 C ) + (7/2) cos2(0 - 9 e - 90°) = (A e - 
7/2) cos2(# — 9 C ). Therefore, while two amplitudes, A c and 7, of ellipticity-only and external shear-only models are 
small, those of realistic ellipticity plus external shear model are large. 

In the case of B2045+265, all models that we considered could not reproduce the observables well (table 4), par- 
ticularly the flux ratios. Some models reduce the x 2 , but we find that they are not physical. However, the system 
shows very large (up to ~ 40%) extrinsic variations of fluxes (Koopmans et al. 2003) on time scales of several months, 
depending on the frequency (Fassnacht et al. 1999). Therefore, it may be premature to conclude that this system 
shows substructure lensing, given that the nature of their variabilities still remains to be understood. One possibility 
is the Galactic scintillation as Koopmans et al. (2003) indicated, and another one is anomalous gravitational lensing 
by the inner structure of the primary lens galaxy. 

4. Conclusion 

We investigated four quadruple gravitational lens system (PG 1115+080, B1422+231, B0712+472, and B2045+265) 
using mass models of SIED or PLF as well as external perturbations described by SIS, or several truncations of the 
expansion of the external perturbation. We performed % 2 fittings with all of the models presented in section 2 in 
order to study the contribution of terms that are usually not included in mass modelings, such as m = 3 multipoles of 
primary lens galaxy and next-order external perturbations. 

We have found that the m = 3 multipoles and/or X3 external perturbation significantly reduce the \ 2 in PG 
1115+080 and B1422+231. In PG 1115+080, not only the external perturbation of the r 3 order (X3), but also an 
internal m = 3 multipole, reduces the \ 2 to statistically acceptable values, but the ellipticity of the galaxy does not. 
Their amplitudes and directions seem to be physical, comparing the results with optical and X-ray observations of 
the galaxy and the group (Impey et al. 1998; Grant et al. 2003); in a (m = 3)+XS model the ratio ^3/00 ~ 0.002 
is small enough to reproduce the observed very spherical galaxy, and in a SIS+(XS+X3) model the ratio (7/7 ~ 0.1 
and the angles 9 1 = —115° and 9 a = —160° can be expected from the positions of the member galaxies and the X-ray 
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brightness contours. A SIS group also provides an acceptable fit. A significant (but small) m — 3 galaxy multipole 
might be induced by interactions with the member galaxies. Regarding B1422+231, several authors concluded that the 
"anomalous" flux ratios reflected the existence of substructure lensing (e.g., Mao, Schneider 1998; Chiba 2002; Bradac 
ct al. 2002; Kochanek, Dalai 2003). However, Keeton, Gaudi, and Petters (2003) found that the cusp relation does not 
reveal a significant anomaly in the system. In fact, our analysis has revealed that a SIED+(XS+X3) model without 
substructure reproduces the fluxes within 5% errors. In addition, the X-ray brightness contours (figure 2 of Grant 
ct al. 2003) imply that the direction of the gradient of the group mass is pa 90° — 180°, and that the shapes of the 
contours are significantly deformed, which are consistent with the results of SIED+(XS+X3) model. In the other 
systems, B071 2+472 and B2045+265, the contribution of the m = 3 multipoles and X3 external perturbations to the 
lens potentials do not significantly improve the model fitting. In B071 2+472, however, we find that the ellipticity and 
shear are well misaligned, 9 C w 56° and 9 1 « 146°, and that the result is very consistent with the lens galaxy and the 
foreground group observed by Fassnacht and Lubin (2002). Therefore, the models including only ellipticity or shear 
may underestimate these (see subsection 3.2). The relatively large shear may arise from the elongated shape of the 
group. 

In summary, the terms of an m = 3 multipole of the primary lens galaxy and of an X3 external perturbation, which 
are often neglected in strong lens analyses, sometimes improve the fit significantly. The contribution of the group 
particularly needs to be studied in more detail. Simple models (in particular one-shear models) may often mislead us. 
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Appendix 1. The Linear Method 



The linear method that we used in this work is as follows. First, we choose a model of Fourier expansion of F(9) 
for a fixed value of /?. Second, for the observed image positions and the relative time delays, we solve a set of linear 
equations [equations (I) and (3)] 

Au = v, (Al) 

where N x M matrix A and N dimensional vector v are functions of observable quantities, such as the image positions 
and the relative time delays, to obtain an initial guess of the parameters ttiincar from which we start to search a 
minimum of the x 2 - Explicitly, the matrix A is expressed as 

/ an a 2 i a 01 a ml (3 ml j n 721 \ 

a il a'21 "ml P'ml 111 721 

. . . . , (A2) 

Axij Ayij 

&0ij (%rnij $m ij "Y 1 ij T2 ij 

v / 

where a lt = cos^, a 2 i = sin 6^, a Qt = fir^ 1 , a mi = fir^ cos m9 t , (5 mi = (3r^ sinmfl,, ju = ncos29i, -f 2 i = ?' 4 sin26' i , 
a'u = — sin#j, a' 2i = cosdi, a' mi = — mrf ~ 1 smm9 i , f3' mi = mr^ 1 cos m6i, j' u = — fjsin2^, 7^ = 7"jcos20j, Ax^j =Xj~Xj, 
Ayij = yi—yj, a 0lJ =rf -r ( -, a mij =rf? cosm^— r^cosmOj, mi j = rf sinm^-r^sinm^-, juj = |(rf cos26>i-r| cos20j) , 
y2ij = ^(rfsia29i-rjaia29j), and T = D~ 1 At lj ([1 + z d ]D d D s /cD ds = Dh~ x ), respectively. The elements of (2i - l)th 
and (2i)th lines are from the lens equations of an image i. The elements of lines that are below (2n)th line arc ones 
of observed time delay for an n imaged system. The vector u is given as 

u T = (x s ,y s ,a ,a m ,b m ,Ji->l2,h). (A3) 

Finally, the vector v is expressed as 

« T = (ri,0,..,i(r?-r?),...). (A4) 

Usually, the number of linear equations N (at most 11 for a quadruple lens system) should be equal or larger than 
that of linear model parameters M. However, we can truncate some of the linear equations so that the number 
of linear equations is the same as that of the parameters, N = M . Thus, the parameters derived from the linear 
method is Mn nC ar = A~ 1 v . We find that this initial guess is close to the actual minimum of the % 2 , provided that the 
observed values have sufficiently small errors. Especially, the Hubble constant, one of the most interesting parameters 
cosmologically, of un ne3X coincides with that of Mbcst within la uncertainties. Figure 1 shows the result of the linear 
method in PG 1115+080, the normalized deviations of ao and h of iiiinear from those of tibest within lcr errors as a 
function of j3 (0.01 for ao and 0.1 for h). Therefore, this method is useful enough to estimate the Hubble constant 
accurately and quickly without solving nonlinear equations. 
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Fig. 1. The normalized deviations of ao and h of ixii near from those of itb cs t m the case of PG 1115+080. We use the time delays 
derived by Barkana (1997). The errors of ao and h of itb cs t are 0-01 and 0.1, respectively. The isothermal model corresponds to 
= 1. 

Our linear method is complementary to that of Evans and Witt (2003) in the sense that our method is applicable 
to arbitrary values of j3 and can derive the external shear as well (the essence of theirs is the fact that for j3 = 1 the 
magnification factor is expressed as a linear equation with respect to the model parameter, 



/i(z)= (1-k) j - 



l7| 5 



An 

— ?Pc 



(A5) 



= { 1 - 2k(1 - 7i cos 26 - 72 sin 20) - | 7 | 2 } , 
when the external shear is assumed). The two methods can be combined in the case of (3 



(A6) 

1; first one estimates 

the external shear with the former and then applies the latter to the considered system in order to reproduce the 
flux ratios without assuming the external shear. However, this combined method sometimes suffers from results of 
unphysical galaxy shapes, such as cross-like ones. 
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